﻿using System;
using System.Diagnostics;
using dotnetcdf;

namespace Faml.Optimization
{
    [OptimizationContainer]
    public static class NetCDFOptimizations
    {
        #region 1d

        [OptimizationForPair("subarray", "netcdf", 2, 0)]
        public static double[] SubarrayNetCDF(string url, string var, int[] indices)
        {
            return SubarrayNetCDF<double>(url, var, indices);
        }

        [OptimizationForPair("subarray", "netcdf_i", 2, 0)]
        public static int[] SubarrayNetCDF_i(string url, string var, int[] indices)
        {
            return SubarrayNetCDF<int>(url, var, indices);
        }

        [OptimizationForPair("subarray", "netcdf_s", 2, 0)]
        public static short[] SubarrayNetCDF_s(string url, string var, int[] indices)
        {
            return SubarrayNetCDF<short>(url, var, indices);
        }

        [OptimizationForPair("subarray", "netcdf", 2, 0)]
        public static double SubarrayNetCDF(string url, string var, int ind)
        {
            return SubarrayNetCDF<double>(url, var, ind);
        }

        [OptimizationForPair("subarray", "netcdf_i", 2, 0)]
        public static int SubarrayNetCDF_i(string url, string var, int ind)
        {
            return SubarrayNetCDF<int>(url, var, ind);
        }

        [OptimizationForPair("subarray", "netcdf_s", 2, 0)]
        public static short SubarrayNetCDF_s(string url, string var, int ind)
        {
            return SubarrayNetCDF<short>(url, var, ind);
        }

        public static T SubarrayNetCDF<T>(string url, string var, int index)
        {
            return SubarrayNetCDF<T>(url, var, new int[] { index })[0];
        }

        public static T[] SubarrayNetCDF<T>(string url, string var, int[] indices)
        {
            int ncid = 0, varid = 0;
            int res, ndimsp = 0;

            res = NetCDF.nc_open(ref url, (int)NetCDF.cmode.NC_NOWRITE, ref ncid);
            if (res != 0)
                throw new Exception("Cannot open NetCDF file");

            try
            {
                int min, max, step;
                bool ordered;
                IndicesUtils.Describe(indices, out min, out max, out step, out ordered);

                res = NetCDF.nc_inq_varid(ncid, ref var, ref varid);
                if (res != 0)
                    throw new Exception("Cannot read NetCDF file");

                res = NetCDF.nc_inq_varndims(ncid, varid, ref ndimsp);
                if (res != 0)
                    throw new Exception("Cannot read NetCDF file");

                if (ndimsp != 1)
                    throw new Exception("Expected for one-dimensional array");

                T[] array = new T[indices.Length];

                if (step > 0)
                {
                    if (typeof(T) == typeof(int))
                        res = NetCDF.nc_get_vars_int(
                            ncid, varid, new int[] { min }, new int[] { indices.Length },
                            new int[] { step }, (int[])((Array)array));
                    else if (typeof(T) == typeof(double))
                        res = NetCDF.nc_get_vars_double(
                            ncid, varid, new int[] { min }, new int[] { indices.Length },
                            new int[] { step }, (double[])((Array)array));
                    else if (typeof(T) == typeof(short))
                        res = NetCDF.nc_get_vars_short(
                            ncid, varid, new int[] { min }, new int[] { indices.Length },
                            new int[] { step }, (short[])((Array)array));

                    if (res != 0)
                        throw new Exception("NetCDF can't read requested variable with a stride function");
                }
                else
                {
                    T[] src = new T[max - min + 1];

                    if (typeof(T) == typeof(int))
                        res = NetCDF.nc_get_vara_int(
                            ncid, varid, new int[] { min }, new int[] { src.Length },
                            (int[])((Array)array));
                    else if (typeof(T) == typeof(double))
                        res = NetCDF.nc_get_vara_double(
                            ncid, varid, new int[] { min }, new int[] { src.Length },
                            (double[])((Array)array));
                    else if (typeof(T) == typeof(short))
                        res = NetCDF.nc_get_vara_short(
                            ncid, varid, new int[] { min }, new int[] { src.Length },
                            (short[])((Array)array));

                    if (res != 0)
                        throw new Exception("NetCDF can't read requested variable in the required region");

                    for (int i = 0; i < indices.Length; i++)
                    {
                        array[i] = src[indices[i] - min];
                    }
                }

                return array;
            }
            finally
            {
                NetCDF.nc_close(ncid);
            }
        }

        [OptimizationForPair("length", "netcdf", 1, 0)]
        public static int lengthNetCDF_d(string url, string var)
        {
            return lengthNetCDF<double>(url, var);
        }

        [OptimizationForPair("length", "netcdf_i", 1, 0)]
        public static int lengthNetCDF_i(string url, string var)
        {
            return lengthNetCDF<int>(url, var);
        }

        public static int lengthNetCDF<T>(string url, string var)
        {
            int ncid = 0, varid = 0;
            int res, ndimsp = 0;
            int dimlen = 0;

            res = NetCDF.nc_open(ref url, (int)NetCDF.cmode.NC_NOWRITE, ref ncid);
            if (res != 0)
                throw new Exception("Cannot open NetCDF file");

            try
            {
                res = NetCDF.nc_inq_varid(ncid, ref var, ref varid);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                res = NetCDF.nc_inq_varndims(ncid, varid, ref ndimsp);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                if (ndimsp != 1)
                    throw new Exception("Expected for one-dimensional array");

                int[] dimids = new int[1];
                int natts = 0;
                res = NetCDF.nc_inq_vardimid(ncid, varid, dimids, ref natts);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                NetCDF.nc_inq_dimlen(ncid, dimids[0], ref dimlen);
            }
            finally
            {
                NetCDF.nc_close(ncid);
            }

            return dimlen;
        }

        #endregion 

        #region 2d
        [OptimizationForPair("length", "netcdf2d", 2, 0)]
        public static int lengthNetCDF2d(string url, string var, int index)
        {
            return lengthNetCDF2d<double>(url, var, index);
        }

        [OptimizationForPair("length", "netcdf2d_i", 2, 0)]
        public static int lengthNetCDF2d_i(string url, string var, int index)
        {
            return lengthNetCDF2d<int>(url, var, index);
        }

        public static int lengthNetCDF2d<T>(string url, string var, int index)
        {
            int ncid = 0, varid = 0;
            int res, ndimsp = 0;

            try
            {
                res = NetCDF.nc_open(ref url, (int)NetCDF.cmode.NC_NOWRITE, ref ncid);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                res = NetCDF.nc_inq_varid(ncid, ref var, ref varid);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                res = NetCDF.nc_inq_varndims(ncid, varid, ref ndimsp);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                if (ndimsp != 2)
                    throw new Exception("Expected for two-dimensional array");

                int[] dimids = new int[2];
                int natts = 0;
                res = NetCDF.nc_inq_vardimid(ncid, varid, dimids, ref natts);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                int dimlen = -1;
                if (index == 0)
                {
                    NetCDF.nc_inq_dimlen(ncid, dimids[0], ref dimlen);
                }
                else
                {
                    NetCDF.nc_inq_dimlen(ncid, dimids[1], ref dimlen);
                }
                return dimlen;
            }
            finally
            {
                NetCDF.nc_close(ncid);
            }
        }

        #endregion

        #region 3d

        #region subarray

        #region for Int32

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[,] SubarrayNetCDF_3d_i(string url, string var, int[] x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[,] SubarrayNetCDF_3d_i(string url, string var, int[] x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[] SubarrayNetCDF_3d_i(string url, string var, int[] x, int y, int z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[,] SubarrayNetCDF_3d_i(string url, string var, int x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[] SubarrayNetCDF_3d_i(string url, string var, int x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[] SubarrayNetCDF_3d_i(string url, string var, int x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int SubarrayNetCDF_3d_i(string url, string var, int x, int y, int z)
        {
           return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_i", 4, 0)]
        public static int[,,] SubarrayNetCDF_3d_i(string url, string var, int[] x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<int>(url, var, x, y, z);
        }

        #endregion

        #region for Double

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[,] SubarrayNetCDF_3d(string url, string var, int[] x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[,] SubarrayNetCDF_3d(string url, string var, int[] x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[] SubarrayNetCDF_3d(string url, string var, int[] x, int y, int z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[,] SubarrayNetCDF_3d(string url, string var, int x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[] SubarrayNetCDF_3d(string url, string var, int x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[] SubarrayNetCDF_3d(string url, string var, int x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double SubarrayNetCDF_3d(string url, string var, int x, int y, int z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d", 4, 0)]
        public static double[, ,] SubarrayNetCDF_3d(string url, string var, int[] x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<double>(url, var, x, y, z);
        }

        #endregion

        #region for Int16

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[,] SubarrayNetCDF_3d_s(string url, string var, int[] x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[,] SubarrayNetCDF_3d_s(string url, string var, int[] x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[] SubarrayNetCDF_3d_s(string url, string var, int[] x, int y, int z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[,] SubarrayNetCDF_3d_s(string url, string var, int x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[] SubarrayNetCDF_3d_s(string url, string var, int x, int[] y, int z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[] SubarrayNetCDF_3d_s(string url, string var, int x, int y, int[] z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short SubarrayNetCDF_3d_s(string url, string var, int x, int y, int z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        [OptimizationForPair("subarray", "netcdf3d_s", 4, 0)]
        public static short[, ,] SubarrayNetCDF_3d_s(string url, string var, int[] x, int[] y, int[] z)
        {
            return SubarrayNetCDF_3d<short>(url, var, x, y, z);
        }

        #endregion

        #region Generic

        public static T[,] SubarrayNetCDF_3d<T>(string url, string var, int[] x, int[] y, int z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, x, y, new int[] { z });
            T[,] res = new T[x.Length, y.Length];
            for (int i = 0; i < x.Length; i++)
            {
                for (int j = 0; j < y.Length; j++)
                {
                    res[i, j] = a[i, j, 0];
                }
            }

            return res;
        }

        public static T[,] SubarrayNetCDF_3d<T>(string url, string var, int[] x, int y, int[] z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, x, new int[] { y }, z);
            T[,] res = new T[x.Length, z.Length];
            for (int i = 0; i < x.Length; i++)
            {
                for (int j = 0; j < z.Length; j++)
                {
                    res[i, j] = a[i, 0, j];
                }
            }

            return res;
        }

        public static T[] SubarrayNetCDF_3d<T>(string url, string var, int[] x, int y, int z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, x, new int[] { y }, new int[] { z });
            T[] res = new T[x.Length];
            for (int i = 0; i < x.Length; i++)
            {
                res[i] = a[i, 0, 0];
            }

            return res;
        }

        public static T[,] SubarrayNetCDF_3d<T>(string url, string var, int x, int[] y, int[] z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, new int[] { x }, y, z);
            T[,] res = new T[y.Length, z.Length];
            for (int i = 0; i < y.Length; i++)
            {
                for (int j = 0; j < z.Length; j++)
                {
                    res[i, j] = a[0, i, j];
                }
            }

            return res;
        }

        public static T[] SubarrayNetCDF_3d<T>(string url, string var, int x, int[] y, int z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, new int[] { x }, y, new int[] { z });
            T[] res = new T[y.Length];

            for (int j = 0; j < y.Length; j++)
            {
                res[j] = a[0, j, 0];
            }

            return res;
        }

        public static T[] SubarrayNetCDF_3d<T>(string url, string var, int x, int y, int[] z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, new int[] { x }, new int[] { y }, z);
            T[] res = new T[z.Length];
            for (int j = 0; j < z.Length; j++)
            {
                res[j] = a[0, 0, j];
            }

            return res;
        }

        public static T SubarrayNetCDF_3d<T>(string url, string var, int x, int y, int z)
        {
            T[, ,] a = SubarrayNetCDF_3d<T>(url, var, new int[] { x }, new int[] { y }, new int[] { z });
            T res = a[0, 0, 0];

            return res;
        }

        public static T[, ,] SubarrayNetCDF_3d<T>(string url, string var, int[] x, int[] y, int[] z)
        {
            int ncid = 0, varid = 0;
            int res, ndimsp = 0;

            res = NetCDF.nc_open(ref url, (int)NetCDF.cmode.NC_NOWRITE, ref ncid);
            if (res != 0)
                throw new Exception("Cannot open NetCDF file");

            try
            {
                int minX, maxX, stepX;
                int minY, maxY, stepY;
                int minZ, maxZ, stepZ;
                bool orderedX, orderedY, orderedZ;
                IndicesUtils.Describe(x, out minX, out maxX, out stepX, out orderedX);
                IndicesUtils.Describe(y, out minY, out maxY, out stepY, out orderedY);
                IndicesUtils.Describe(z, out minZ, out maxZ, out stepZ, out orderedZ);

                res = NetCDF.nc_inq_varid(ncid, ref var, ref varid);
                if (res != 0)
                    throw new Exception("Cannot read NetCDF file");

                res = NetCDF.nc_inq_varndims(ncid, varid, ref ndimsp);
                if (res != 0)
                    throw new Exception("Cannot read NetCDF file");

                if (ndimsp != 3)
                    throw new Exception("Expected for three-dimensional array");

                T[, ,] array = new T[x.Length, y.Length, z.Length];

                bool regular = (stepX > 0) && (stepY > 0) && (stepZ > 0);
                if (regular)
                {
                    int xs = x.Length;
                    int ys = y.Length;
                    int zs = z.Length;
                    T[] src = new T[xs * ys * zs];

                    if (typeof(T) == typeof(int))
                        res = NetCDF.nc_get_vars_int(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { x.Length, y.Length, z.Length },
                            new int[] { stepX, stepY, stepZ }, (int[])((Array)src));
                    else if (typeof(T) == typeof(double))
                        res = NetCDF.nc_get_vars_double(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { x.Length, y.Length, z.Length },
                            new int[] { stepX, stepY, stepZ }, (double[])((Array)src));
                    else if (typeof(T) == typeof(short))
                        res = NetCDF.nc_get_vars_short(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { x.Length, y.Length, z.Length },
                            new int[] { stepX, stepY, stepZ }, (short[])((Array)src));
                    
                    if (res != 0)
                        throw new Exception("NetCDF can't read requested variable with a stride function");

                    int index = 0;
                    for (int i = 0; i < xs; i++)
                    {
                        for (int j = 0; j < ys; j++)
                        {
                            for (int k = 0; k < zs; k++)
                            {
                                array[i, j, k] = src[index++];
                            }
                        }
                    }
                }
                else
                {
                    int xs = (maxX - minX + 1);
                    int ys = (maxY - minY + 1);
                    int zs = (maxZ - minZ + 1);
                    T[] src = new T[xs * ys * zs];

                    if (typeof(T) == typeof(int))
                        res = NetCDF.nc_get_vara_int(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { xs, ys, zs }, (int[])((Array)src));
                    else if (typeof(T) == typeof(double))
                        res = NetCDF.nc_get_vara_double(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { xs, ys, zs }, (double[])((Array)src));
                    else if (typeof(T) == typeof(short))
                        res = NetCDF.nc_get_vara_short(ncid, varid, new int[] { minX, minY, minZ },
                            new int[] { xs, ys, zs }, (short[])((Array)src));

                    if (res != 0)
                        throw new Exception("NetCDF can't read requested variable in the required region");

                    int l, m, n;
                    for (int i = 0; i < x.Length; i++)
                    {
                        l = (x[i] - minX) * ys * zs;
                        for (int j = 0; j < y.Length; j++)
                        {
                            m = (y[i] - minY) * zs;
                            for (int k = 0; k < z.Length; k++)
                            {
                                n = z[i] - minZ;
                                array[i, j, k] = src[l + m + n];
                            }
                        }
                    }
                }

                return array;
            }
            finally
            {
                NetCDF.nc_close(ncid);
            }
        }

        #endregion

        #endregion

        #region length
        [OptimizationForPair("length", "netcdf3d", 2, 0)]
        public static int lengthNetCDF3d(string url, string var, int index)
        {
            return lengthNetCDF3d<double>(url, var, index);
        }

        [OptimizationForPair("length", "netcdf3d_i", 2, 0)]
        public static int lengthNetCDF3d_i(string url, string var, int index)
        {
            return lengthNetCDF3d<int>(url, var, index);
        }

        [OptimizationForPair("length", "netcdf3d_s", 2, 0)]
        public static int lengthNetCDF3d_s(string url, string var, int index)
        {
            return lengthNetCDF3d<short>(url, var, index);
        }

        public static int lengthNetCDF3d<T>(string url, string var, int index)
        {
            int ncid = 0, varid = 0;
            int res, ndimsp = 0;

            try
            {
                res = NetCDF.nc_open(ref url, (int)NetCDF.cmode.NC_NOWRITE, ref ncid);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                res = NetCDF.nc_inq_varid(ncid, ref var, ref varid);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                res = NetCDF.nc_inq_varndims(ncid, varid, ref ndimsp);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                if (ndimsp != 3)
                    throw new Exception("Expected for three-dimensional array");

                int[] dimids = new int[3];
                int natts = 0;
                res = NetCDF.nc_inq_vardimid(ncid, varid, dimids, ref natts);
                if (res != 0)
                    throw new Exception("Cannot open NetCDF file");

                int dimlen = -1;
                NetCDF.nc_inq_dimlen(ncid, dimids[index], ref dimlen);
                return dimlen;
            }
            finally
            {
                NetCDF.nc_close(ncid);
            }
        }

        #endregion 

        #endregion
    }
}
